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Abstract 

Multilevel (ML) algorithms for eigenvalue problems are often faced with several types 
of difficulties such as: the mixing of approximated eigenvectors by the solution process, the 
approximation of incomplete clusters of eigenvectors, the poor representation of solution on 
coarse levels and the existence of close or equal eigenvalues. Algorithms that do not treat 
appropriately these difficulties usually fail, or their performance degrades when facing them. 
These issues motivated the development of a robust adaptive ML algorithm which treats 
these difficulties, for the calculation of a few eigenvectors and their corresponding eigenval- 
ues, presented in this paper. The main techniques used in the new algorithm include: the 
adaptive completion and separation of the relevant clusters on different levels, the simul- 
taneous treatment of solutions within each cluster, and the robustness tests which monitor 
the algorithm’s efficiency and convergence. The eigenvectors’ separation efficiency is based 
on a new ML projection technique generalizing the Rayleigh Ritz projection, combined with 
a novel technique, the backrotations. These separation techniques, when combined with an 
FMG formulation, in many cases lead to algorithms of 0{qN) complexity, for q eigenvectors 
of size N on the finest level. Previously developed ML algorithms are less focused on the 
mentioned difficulties. Moreover, algorithms which employ fine level separation techniques 
are of 0{q^N) complexity and usuaUy do not overcome all these difficulties. Computational 
examples are presented where Schrodinger type eigenvalue problems in 2-D and 3-D, hav- 
ing equal and closely clustered eigenvalues, are solved with the efficiency of the Poisson 
multigrid solver. A second order approximation is obtained in 0{qN) work, where the total 
computational work is equivalent to only a few fine level relaxations per eigenvector. 


*This research was made possible in part by funds granted to the second author through a 
fellowship program sponsored by the Charles H. Revson Foundation and in part by the National 
Aeronautics and Space Administration under NASA Contract No. NASl-19480 and NASl-18605 
while the author was in residence at ICASE, NASA Langley Research Center, Hampton, Va 23681 
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1 Introduction 


Large scale eigenvalue problems (EP) arising from physics, chemistry and engineering 
often have special features which are not always exploited by eigenvalue solvers, such 
as: the EP can be approximated on several levels; only few eigenvalues and eigenvec- 
tors are sought. These features can be exploited naturally by multilevel (ML) solvers 
as has been done successfully by several authors (see for example [3] and [12]). Such 
solvers generally involve the discretization of the problem on a sequence of levels, 
relaxations employed on all levels, transfers of solutions and residuals from fine to 
coarse levels, interpolation of correction from coarse to fine levels, and often a fine 
level eigenvector separation technique. 

Multilevel methods for EP sometimes encounter difficulties which make their ro- 
bustness and efficiency questionable. Often such difficulties are due to several reasons 
which we classify as: (i) mixing of eigenvectors by the used procedures, (ii) incom- 
pleteness of a treated cluster, (iii) incompatible fine-coarse level representation of a 
cluster. Procedures such as relaxations, transfers and ML cycles can introduce or am- 
plify an eigenvector in the error component of an approximated eigenvector. This we 
refer as eigenvector mixing. For example, it is known that ML procedures mix Fourier 
components which are the eigenvectors of many discretized differential operators. If 
the eigenvectors of an iteration operator do not coincide with the eigenvectors of the 
problem to be solved, then the iteration may mix the problem’s eigenvectors. This dif- 
ficulty occurs especially when clusters of eigenvectors with close or equal eigenvalues 
are approximated. Usually it is treated by simultaneous separation techniques, e.g., 
by a Rayleigh Ritz type projection. However, if not all eigenvectors which are mixed 
by the involved procedures are approximated, e.g. the clusters are not complete, then 
one can expect the separations to be inaccurate and inefficient. Thus a cluster will 
be called complete relative to a procedure if it contains a whole set of eigenvectors 
which are mixed by the procedure. Difficulties related to incomplete clusters may be 
treated by completing the clusters and processing them simultaneously. The incom- 
patibilty of clusters representation on different levels is another difficulty which has to 
be identified and taken care of. In particular, not all clusters can be approximated on 
an arbitrarily coarse level and the eigenvalues clustering may differ on various levels. 
These difficulties and their remedies suggest the following conclusions: (i) Clusters of 
eigenvectors should be treated simultaneously using a separation technique like the 
Rayleigh-Ritz projection, (ii) Clusters have to be completed, (III) Different clusters 
should be differently treated on different levels. This suggests that previously de- 
veloped algorithms may fail in some standard situations. These include nonadaptive 
algorithms, algorithms which treat the eigenvectors sequencially in clusters, which do 
not complete the relevant clusters, or which do not take into account the inter-level 
eigenvector mixing. Even in cases when such algorithms work, their efficiency may be 
improved by an adaptive treatment or by coarse level separation techniques as used 
in this work. 
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This paper focuses on a more robust and efficient algorithm for the calculation of 
a few eigenvalues and their corresponding eigenvectors. Its development was guided 
by the above mentioned difficulties and their remedies. Beside its robustness, the 
algorithm achieves a better computational complexity than previously known ML 
eigenvalue algorithms which use fine level projections. The robustness of the present 
approach is based on the adaptive completion and separation of the relevant clusters 
on different levels; the simultaneous treatment of solutions for each cluster; and ro- 
bustness tests which monitor the algorithm’s convergence and efficiency. A central 
efficiency feature of the algorithm presented here results from the newly developed 
Generalized Rayleigh Ritz (GRR) projection and backrotation (BR), which employ 
the projection on coarse levels, adaptively with the cluster involved. This reduces in 
many Ceuses the most time consuming part of the algorithms, namely, the 0{q^N) fine 
level separation work, to 0{qN) for q eigenvectors of size N on finest level. 

These ideas are combined in an FMG algorithm which first solves the problem 
on coarse levels, then interpolates the solutions to finer levels where they serve as 
an initial approximation to the corresponding finer level problems. On the currently 
finest level, clusters are identified and tested for completeness, completed if necessary, 
and improved by ML cycles using coarser levels. The eigenvalue equations is relaxed 
on each level followed by FAS transfers. Generalized Rayleigh Ritz projection and 
backrotations (GRR-BR) are employed on coarse levels usually, in order to separate 
eigenvectors within their clusters, and to keep the coarse level representation of the 
solutions as close as possible to the fine level solutions. This is done adaptively for 
different clusters on appropriate lev^._On the level on which the algorithm starts, 
only a part of the sought eigenvectors are approximated usually, and more eigenvectors 
with their corresponding eigenvalues are added on finer levels. 

In the examples presented here for the Schrodinger eigenvalue problem in 2-D 
and 3-D, two to four fine level relaxations per eigenvector were performed. Equal 
eigenvalues were calculated with more than ten decimal places and accurate results 
were obtained for very close eigenvalues as well. A second order approximation is 
obtained by 1-FMG-V(1,1) in OlqN) work, for q eigen vec tors of size N on the finest 
level. V -r 

The present approach can be extended to nonlinear eigenvalue problems, an ex- 
ample being presented in Costiner and Ta’asan [6]. 

We refer to the early works of Hackbusch [9] McCormick [17] , Bank [1] for theory 
on multilevel eigenvalue solvers and first algorithms. A sequencial ML algorithm for 
linear eigenvalue problems performing the projection on fine levels is presented in 
Brandt McCormick and Ruge [3]. More theory and algorithms on ML EP may be 
found in Hackbusch [12]. Zaslavski in [23], [24], uses an adaptive algebraic correction 
scheme cycle to compute the first eigenvector and its eigenvalue for the multigroup 
neutron diffusion equation. The elements of such a ML cycle, modified to the FAS 
form, can be used in the algorithm presented here. Our approach differs from previous 
ML approaches mainly by the following issues: emphasis on robustness, simultaneous 
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cluster processing, cluster completion, ML separation by GRR and BR, treatment of 
eigenvector mixing, treatment of close and equal eigenvalues. The ML projection idea 
was first introduced by Ta’asan [21]; backrotations were introduced in Costiner [5]. 
The combination of our technique with domain decomposition techniques is natural 
but was not analysed yet. For some domain decomposition techniques for eigenvalue 
problems see for example Bourquin and d’Hennezel [2], and Luo [14]. A review article 
on single level large-scale complex eigenvalue problems, containing many references, is 
Kerner [13]. For a theory on Ritz projections and on algebraic eigenvalue problems 
we refer to Parlett [19], Wilkinson [22], Golub and Van Loan [8]. The single 
level technique to obtain the eigenvectors by relaxations and projections is refered in 
different places as subspace, simultaneous or Ritz iterations. We refer to Nikolai [18], 
Rutishauser [20] and McCormick [17] for a single level algorithm and mathematical 
foundations. 

The paper is organized as follows. Section 2 presents the GRR ML projection, the 
backrotations and the multilevel cycle. Section 3 presents and discusses the adaptive 
techniques such as the robustness tests, the cluster completion, and the adaptive 
FMG. Section 4 presents computational examples. 


2 Multilevel Projection Techniques 

One of the key elements of our algorithm is the ML projection used to separate 
eigenvectors corresponding to closely clustered and equal eigenvalues. To motivate 
the new projection method consider a fine level problem 

- (7,A = 0 (2.1) 

where A = diag{\i, . . . , Xg) contains on diagonal the q sought eigenvalues corre- 
sponding to the desired eigenvectors which are the columns of Uh- Assume that Uh. 
consists of linearly independent combinations of eigenvectors belonging to the sought 
subspace. In this case, a Rayleigh-Ritz (RR) projection provides a A and & q x q 
invertible matrix E such that A and Uh = UhE are solutions of (2.1). The FAS coarse 
level equation after performing the projection on the fine level becomes 

A2hU2hE - UihEA = T^hE (2-2) 

where = A 2 hIh'’Eh ~ Il^AhUh is the usual FAS right hand side. The form of 
the coarse level solution U 2 h,E = l^^UhE suggests that one may obtain E and A on 
the coarse level from 2.2. This can be done by a generalized version of the (RR) 
projection, presented next. 
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2.1 Generalized Rayleigh Ritz Projections and Backrota- 
tions 

Solutions [E, A) of (2.2) may not exist if (//, is not a basis of span[Uh)- However, 
as in the usual (RR) projection, if approximates a basis of span{Uh)i one finds 
E and A such that the projection of the residual of (2.2) on span{U 2 h) is minimized, 
i.e., solving a. q x q generalized eigenvalue problem as: 

iU2kf{A2kU2h - T^h)E = {U2kfU2hEA (2.3) 

The process of solving for (£, A) given (A 2 A, Uih.T^h) is denoted by 

[E, A) ^ GRR(A2a, U2h, T4) (2.4) 

The above projection is. refered later on as generalized Rayleigh Ritz Projection 
(GRR) or as multilevel ML projection when several levels are involved. The GRR 
projection cannot be directly combined with the usual FAS correction Uk = Uh 

since it will change an exact fine level solution f//,, e.g. if E is not 
the identity but rather a permutation matrix. This difficulty can be solved by using 
a modified FAS correction such as 

f/A = + 4(f/2A£^ - I'fHlhE) (2.5) 

Note that ( 2.5 ) would lead to 0{q^N) operations, equivalent to a fine level projection 
work. Thus it is desirable to replace ( 2.5 ) with more efficient techniques. Other 
difficulties may occur for degenerate subspaces^ien any matrix can serve cis a solution 
for E, thus, mixing and destroying orthogonality of fine level solutions for example. 

A natural technique to fix these difficulties is introduced and used in combination 
with the multilevel projection. It is refered to as backrotation because of its geo- 
metrical meaning. As suggested by the above discussion, for degenerate subspaces 
the backrotation should produce block diagonal submatrices in E which are close 
to the identity matrix of the appropriate dimension and eliminate permutations of 
eigenvectors. A backrotation step will be further denoted by 

(£, A) BackRotation(£^, A) (2-6) 

A particular backrotation algorithm is: 

BackRotat ion 
Input {E,A) 

1 ) Sort the eigenvalues of A and 

permute the columns of E accordingly 

2) Determine the clusters of eigenvalues of A 

to be considered degenerate, and 

determine the clusters to be considered nondegenerate 
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3) For each diagonal block in E 

associated with a nondegenerate cluster do: 

Bring to the diagonal the dominant elements of the block by 

permuting the columns of E, 

and correspondingly the diagonal of A. 

4) Let F be a block diagonal matrix 

whose diagonal blocks are the diagonal blocks of E, 
corresponding to the clusters. 

Replace each diagonal block which does not correspond 
to a degenerate cluster by the corresponding identity matrix 

5) Set E = EF~K 

6) Change the signs of columns of E 

to get positive elements on diagonal. 

7) Normalize the columns of E. 

Output {E, A) 


2.2 Multilevel Combined Cycles 

In general, not all required eigenvectors can be well approximated on an arbitrarily 
coarse level. For example, fine level discretized Laplace operator eigenvectors corre- 
sponding to larger eigenvalues may not be representable on coarse levels. Moreover, 
the eigenvector cluster structure may differ on various levels. However, the efficient 
solution of the eigenvectors corresponding to small eigenvalues may often be done 
using coarse levels. A major difficulty is that solutions belonging to a cluster are 
often mixed by procedures, i.e. the procedures may regard linear combinations of 
solutions as a solution. This obstacle can be overcomed by simultaneously treating 
all the cluster’s eigenvectors and separating them on different levels. 

Efficiency and convergence considerations require that the GRR projection should 
be done for different clusters on different levels. Moreover, the coarsest level used to 
treat a given cluster may not coincide with the level on which the GRR projection 
is done. Thus, the full algorithm depends on some parameters associated with each 
cluster, that determine the flow of the algorithm for that cluster. These parame- 
ters are determined adaptively during the solution process (as explained in the next 
section). 

Following is a description of a basic ML cycle used in the adaptive algorithm 
presented in the next section. 

Let q’ eigenvectors be approximated by j clusters on level k: 

U, = (Ul,...,Ui) (2-1) 

where, each Ul, approximates (/^, the solution of 

A,Ui = UiA<+n i = (2.2) 
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where each Ul is a matrix whose columns are eigenvectors and A' is a diagonal ma- 
trix whose diagonal elements are the corresponding eigenvalues. Usually, on the finest 
level, k = m, Tk = {Tl,... ,Tl) = (0, . . . ,0). Denote A = diag{A\. . . , A-'). For each 
cluster Ul let be the level on which the GRR-BR projection is done, and the 
coarsest level used in the ML process for this cluster. Here it is assumed that 
Denote Ip = (/’, . . . , l^), L = (/],.. . Ji) and by II a function transfer from level k to 
level j. For improving a given approximation {Umi A,T,„), a multilevel cycle consist- 
ing in a sequence of cycles for each cluster in turn, is: 

{Umi ^1 Tmi ^pi ^C1 9 O * — CL-MLP (tTI, Ami ^^mi A, Tmj ^pi ^ci 9 ) 

For i = 1, . . . , j do: ' 

For L- = rn, ..., /^ do: 

Repeat Nl Times: 

If li, = k then {Ui, A\Tl) 4- GRR-BR(m, Aki Ul A', T/, k) 

Ui <-ReIax {m, AkiUl, A\Tl,k,ll) 

If k > li then: 

Set k — k — li 

ui = it+^uu„ 

Ti = IL, {TU, - Am U’m ) + AtUi 

End 

For L- = /*,... , rn do: - : . 

If (k > /<) Correct Ui = (/j + ILM., - it'Ul) 

Repeat Nl Times 

Ui <r-Relsoc {m,AkiUl,A\Tl,k,ll) 

If = jfc then {Ui, A\ Ti) e- GRR-BR(m, Aki Ui, N,Tl, l;, k) 

End 

End 


The GRR-BR algorithm used above is the following: 

(^/’ , A^ ^*)>- GRR-BR{m,A,,U^ 

Perform - _ 

{E,A')^G^{Ak,Ui,Tl^ 

(F,A’) <— BackRotatation(F, A') 

ui = uiE 

n = TiE 


Observe that in CL-MLP the clusters are treated sequentially and within each 
cluster the solutions are treated simultaneously in the ML cycle. A simultaneous 
cycle for several clusters is obtained by grouping the clusters into a single larger 


6 



cluster and applying CL-MLP to it. This may be used to improve separation between 
clusters. Observe that if for each cluster the GRR-BR projection is performed on 
the finest level, the algorithm still requires less work than an algorithm performing 
the fine level projection for all clusters simultaneously. Moreover, if mixings occur 
on coarse levels, one may expect an algorithm using fine level separation to have 
a poor efficiency; a coarse level separation usually improves the efficiency in such 
cases. For well separated eigenvalues the projection may not be needed except at 
initial coarse level stages of the FMG, later the eigenvalues determine the separation 
of eigenvectors via the multilevel cycles. The same holds for well separated clusters 
which do not need a simultaneous separation. This is especially useful for a larger 
number of eigenvectors, belonging to well separated clusters (e.g. already for 10 
eigenvectors the improvement can be noticeable). Accurate normalization, if needed, 
can be performed as the last step only on the finest level. 

Several parameters in this algorithm have to be specified, such as coarsest levels 
Ic and Ip, and numbers of iterations. These parameters are chosen by the adaptive 
algorithm discussed in next section. The choice of Ip and Ic depends on cluster and 
stage of the algorithm. The number of iterations Nk can be deduced for simple cases 
based on Fourier or 2-level cycle analysis or in general cases, by robustness tests 
performed during execution. The number of relaxations can vary with level. In the 
computational tests one or two relaxations per fine level passing were performed. 

sought 

3 Adaptive Multilevel Algorithms 

For robustness, the construction of an adaptive version of the CL-MLP is essential 
since the grouping of eigenvectors into clusters, the sizes of the clusters, the coarsest 
level corresponding to a given cluster are not known in advance, usually. Eigenvectors 
belonging to clusters usually are mixed by different procedures. This deteriorates the 
algorithms’ efficiency and often prevents convergence. A typical difficulty occurs 
when a procedure approximates only several eigenvectors of a cluster. Then the 
nonapproximated eigenvectors usually are the dominant components of the errors 
which are hard to eliminate due to mixing. This suggests to complete the clusters and 
separate the solutions within clusters whenever necessary. Simultaneous techniques 
treating at a time all solutions belonging to a complete cluster, can be easier coupled 
with separation techniques at different stages, thus acquiring better efficiency than 
sequencial techniques which hardly avoid difficulties due to mixing. 

In the adaptive algorithm, the clusters are tested for completion and completed. 
The cluster completion is tested on all current finest levels and performed on several 
levels since the structure of clusters can differ on different levels. 

The full multilevel solver described below starts on coarse levels, and solves there 
for as many eigenvectors as possible. Then it uses those as an initial approximation 
for finer level solutions where more eigenvectors are added if needed. 
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Two important parts of this algorithm are the completion and addition of clusters. 
The completion of a cluster is done by adding in turn a new vector and improving it 
by multilevel cycles. An approximate eigenvalue is associated with this eigenvector, 
by a Rayleigh quotient. If the eigenvalue is close to the cluster then the new vector 
is added to the cluster. If the found eigenvalue does not belong to the cluster then 
the cluster is considered complete. The convergence of the additional eigenvector is 
not sought. At the end, the complete cluster is improved by several CL-MLP cycles. 

The addition of a new clusters is usually done in the first stages of the algorithm 
when not enough complete clusters are found. This is performed by the cluster 
completion algorithm described next. 

Denote by dj the current dimension of the cluster (/^. The cluster completion and 
cluster addition algorithms are given by: 


{j, Ui,A\ TiX, h, q') ^Complete-Cluster(i, Ul A\ T/, q') 

Until (Cluster-Completion-Test = TRUE) Do 
Choose random <j> ^ 

Until < Ak<f>, (f) > f < (j>,(f>> and residuals stabilize Do: 

((S, a;,„, n,0, 1) ^cl-mlp(J:,.^, Aj.„, n', 

Separate ^ from (f/fc, ..., (/^) 

Set A =< Ak(j>, (f> > j < <f>,4>> — - — 

Vl-(ul4)_ 

A-' diag{\^ , A) 

q' = q' dj = dj + 1 

End : ■ 

Perform {Ul, A^TU^Jl^dj) ^CL-MLP{k,Ul,A\Tl,P^,ll,dj) 


{j, Uk, A, Tk, Ip, Ic, q') ^ Add-Cluster(j, Uk, A, Tk, Ip, h, q') 

Set j = j + 1 

(i, Ul, A\ Tl, Ip, h, q') e-Complete-Cluster(j, Ul, AU Tl, V^, l^, q') 
Set Uk = {Ul,...,Ui), A = (A^...,A^■) 


The Separation of <j> from the other eigenvectors may be performed by orthogo- 
nality, projection or by a ML cycle. 

Another crucial part of the algorithm are the robustness-tests. Robustness tests 
are techniques which find parameters to be used in a certain procedure for given 
data such that the procedure will be convergent and will be a.s efficient as possible. 
Typically, the procedure is tested over a set of data, or information from intermediate 
results is used in order to reduce the work involved in testing. Next, for simplicity. 
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we consider only the robustness test which provides the parameters {lp,lc) for the 
CL-MLP cycle. This ensures the proper values for the parameters Ic, Ip. in order 
to get convergence and the best attainable efficiency for CL-MLP. This is done by 
examining for clusters the different possible values for Ic, Ip. The use of an FMG 
algorithm suggests that these parameters for a given cluster will stabilize with the 
refinement process. Thus, a notion of a stabilized cluster is introduced in order to 
save unnecessary work in looking for new values for Ic, Ip for stabilized clusters. 

A complete cluster on level L is called stabilized if it corresponds to a complete 
cluster from level L— 1 or L + 1 in the sense of the number of eigenvectors in the cluster, 
the values of the eigenvalues and the eigenvectors approximation. For stabilized 
cluster, corresponding to a coarser level stabilized cluster, we take the Ic, Ip values 
from the corresponding coarser level cluster. For not stabilized clusters, which would 
exist usually on coarse levels only, we perform a search to obtain best values for Ic, Ip. 
This is done by performing ML cycles with different choices of these parameters, and 
choosing the ones that perform best. This is summarized in the next algorithm. Such 
tests are inexpensive when performed on coarse levels, and often lead to significant 
fine level work savings. Moreover, such tests are essential to ensure the algorithm’s 
convergence. 

Convergence is always obtained since at least the single level cycle converges, being 
a subspace iteration algorithm [8]. 

Denote by (/^, l^, m) the Ip and Iq parameters, for the level m cycle for the cluster 
j, and by := n{C L - M LP{m, A,n,Um, ^,Tm,ll,li,/)) the convergence 

rate (measured by the residual decrease) of the CL-MLP cycle for cluster j on level 
m, using the parameters {Ip, Ic). The following algorithm updates {lp,lc) on level m: 

{lp,lc,m) <- Robustness- Test {m, A,n,U,n, A,T,n,lp,lc,q') 

For j Do: 

If - aL_.,| < c ) 

then 

- 1 ) 

else 

If (|a;„ — I > e ) or if (A^ ) is not approximated 
then 

Solve for 

min,j fi{lp, lc,m) : Vq < Ip < m, 
else 

= {Pp,li,m-l) 

endif 

endif 

Observe that since the initial values for li, on any level are taken from the next 
coarser level, the search needed in the Robustness- Test is over just a few choices of 
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the paranieters. 

The full ML algorithm uses as building blocks the CL-MLP, Add-Cluster, Complete- 
Cluster and Robustness-Test algorithms described before. It is defined as: 

Ad[aptive-FMG(m, A) 

Set k = 1, q' = 0, j = 0, Pj, = A;, = k 

Until {q' > q or q' > adim/t) Perform 

(i, Uk, A, Tk, Ip, Ic, q') <-Add-Cluster{j, Uk,A,Tk, Ip, h, q') 

{Uk, A,Tk,lp,lc,q') ^C\.-MLV{k,Uk,A,Tk,lp,lc,q') 

Until k > m Do: 

If k < m then: 

SetA: = A:-fl, Uk = Ik-ilh-i, Ifc = 0 
endif 

If (q' > (l) then: 

if (Cluster-Completion- T'est=TllUE then: 

{Uk, A, Tk, lp,lc, q') ^CL-MLP(A;, Uk, A, Tk, Ip, h, q') 

Else ' 

{j, ui, A\ TI, Pp, P^, g') ^Complete-Cluster(i, Ul, A\Tl , P^, P^, q') 
{Uk,A,Tk,lp,lc,q')^CL-MLP{k,Uk,\,Tk,lp,lc,q') 
endif 

Else-"-. 

Until {q' > q or q' > adim*) Pecfprtn 

{j, Uk, A, Tk, Ip, Ic, q') <-Add-Cluster(j, Uk, A, Tk, Ip, h, q') 
[Uk,A,Tk,lp,h,q') ^Cl.-MLP{k,Uk,\,Tk,lp,lc,q') 

endif 

3.1 Storage and Complexity 

For the Adaptive-FMG algorithm, storage is required for the q eigenvectors of size N 
on the finest level, the potential and the corresponding right hand sides, on all levels, 
giving an overall estimate of memory of order 0(3(Af -b 1)) for problems in 2-D and 
3-D. The FMG work requires 0{N) operations per eigenvector. The work performed 
on coarsest levels should be added to these estimates. In the performed tests where 
a few eigenvectors were sought, the coarse level work was usually a fraction of the 
finest level work. If accurately zero scalar products are needed on finest levels then 
orthonormalizations or projections may be required within the finest level degenerate 
or close clustered eigenspaces. However, as can be seen in the computational ex- 
amples, accurate orthogonality inside degenerate clusters may be obtained by coarse 
level separation also. 

In the computational examples presented here, a complexity of 0{qN) is obtained. 
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4 Computational Examples and Discussions 

We have chosen several examples in which we focus on different aspects of the method. 
The first example is of an adaptive algorithm in which we take a case of special 
difficulties in terms of clustering of eigenvalues, and correspondence of these clusters 
between levels. The second example shows that it is enough to treat the clusters in a 
sequential manner, and separation has to be done within each cluster only to obtain 
good convergence and accurate separation. 

The third example shows that the new multilevel projection technique may be 
performed on coarsest levels, even in cases with closely clustered or equal eigenvalues, 
thus reducing the computational work significantly. A last example shows that the 
same efficiency is obtained for problems in 3-D as for problem in 2-D. 

In all examples the periodic boundary conditions Schrodinger eigenvalue problem 

(A — V)u = Au (4.1) 

defined on D = [ 0 , 0 ]“^ (d=2 or 3) where a = 27r/10, was considered. The potentials 
are chosen such that the distribution of eigenvalue present special difficulties. A finite 
difference discretization on cartesian grids is used. 

Example I: Adaptive Algorithm 

As an example for the adaptive algorithm we have chosen a potential which deter- 
mines a distribution of eigenvalues which present special difficulties. Not only closely 
clustered eigenvalues are present but also the correspondence of the clusters between 
levels is not appropriate, as explained later. 

The potential for this case is V (x, y) = 5-l-3stn(10x) and the results are presented 
in Tables 1 and 2. The first q = 12 eigenvalues were required, and have been approx- 
imated using an adaptive 1-FMG-V(1,1) algorithm which uses as its coarsest level a 
employs 4x4 coarsest grid. 

The ith eigenvalue and eigenvector will be denoted next by A,- and u,-. The boxes 
in Table 2 show the clusters of close or equal eigenvalues (with (-) sign) found by 
the algorithm (the formats are chosen to outline the equal digits). As can be seen 
there the cluster structure on the different levels is not the same. Particularly, level 2 
cluster structure differs from the level 1 cluster structure. Observe that the cluster of 
6 eigenvalues on level 1 (Ae ~ An), with multiplicities 1 - 4 - 1, has no correspondence 
on level 2. 

For demonstrating the adaptive flow of the algorithm we give a full history of one 
specific run, for the potential discussed here. 

The particular algorithm described in section 3 is used. 

The algorithm started on level 1 adding eigenvectors until the cluster containing 
Ai 2 was completed. Note that Aie, the last eigenvector that was found belongs to the 
next cluster, confirming the completeness of the last sought cluster. Observe that 
on level 1, A 12 belongs to a cluster consisting of two degenerate subspaces, each of 
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dimension 2, and the eigenvalues corresponding to these degenerate subspaces are 
close to within 0(10““*) relative difference. 

The relevant eigenvectors Ui, . . . ,ui5 were interpolated to level 2 where they pro- 
vided initial guesses for the level 2 algorithm. Here the completion of clusters restarted 
but this time working with the cluster structure from level 1 and using two level cycles. 
A test was done for the efficiency of a simultaneous cycle with fine level projection. 
The cycle was performed to provide first approximations of the level 2 eigenvalues. 
The cluster structure and eigenvalues obtained were compared with the ones of level 
1. Since the agreement was not satisfactory, except for Ui, a cluster completion algo- 
rithm started with v^- The completion continued until the complete cluster containing 
the last sought eigenvector was obtained, (e.g. for level 2, the desired belongs to 
the cluster Uio — uis- The completion was ensured by the far value of A14). Then the 
relevant eigenvectors were updated by one or a few cycles. 

The solution obtained oii level 2 was interpolated to level 3 where a cluster comple- 
tion test was satisfied only by the first cluster, Ui. The cluster completion algorithm 
was applied to the remaining eigenvectors (using cycle robustness tests and the cluster 
completion tests). This resulted in few cycles per eigenvector. The parameters Ic and 
Ip were found in the following way: 1) for cluster 1, (ui), the values were obtained 
from previous level since this cluster was stabilized from the beginning; 2) For cluster 
2 and 3 (ue ~ ^^9 a-^d Uio — uia ) ic and Ip were taken from level 2 values since these 
clusters resulted stabilysed after the cluster completion; 3) Robustness tests were 
used for cluster 4 since the eigenvalues Aio — A13 on level L=3 and the next coarse 
level L=2 were not enough close. Then one cycle (here V(l,l)), was performed for 
each cluster. 

The level 4 cluster completion test was satisfied by the first 3 clusters, eigenvalues 
U] — Vg and their parameters were taken from level 3. The cluster completion algo- 
rithm was again applied to cluster 4, (ujo — U13), a few cycles being sufficient, then 
the parameters were taken the corresponding level 3 ones since the cluster resulted 
stabilysed. One cycle was performed for each cluster. 

The level 5 cluster completion test was satisfied by all relevant eigenvectors (ui — 
U13). One cycle (V(l,l)), was performed for each cluster. The Ic and Ip for the separate 
clusters, in the final cycles, on levels 3, 4, 5 were found as: for vj: Ic = ip = 1, for 
the other clusters, containing u-2, . . . ,vj3 ic = ip = 2 were obtained, ( a test for the 
asymptotic convergence rate, for cluster Uio — uis may lead to h = ip = 3, but such a 
test was not used in this run). 

The additional last eigenvector obtained in the cluster completion test, used just 
to ensure that the previous cluster was complete, was not needed and not used in 
further steps. Usually its convergence was poor since the algorithm didn’t separate 
it from the next eigenvectors in its cluster e.g. on level 2, to separate A14 from the 
next 7 eigenvectors with close eigenvalues. 

Observe that a second order approximation and a good convergence rate of order 
O(10~^) for the first cycle are obtained. A simultaneous cycle for all clusters with 
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separation on the coarsest common level for all clusters (here level 2 ) would improve 
the efficiency of first cycle but was not needed. (This also would improve the scalar 
products which resulted of order 10"'* after first FMG cycle, in this case, accurate 
orthogonality being obtained by the algorithm described next). 

Observe that this algorithm is of order 0{qN) if one does not use fine level sep- 
aration inside the clusters. The adaptive coarse level work on levels 1, 2, took ap- 
proximately 1/6 of the computer time and on levels 1, 2, 3, approximately 1/4 of the 
computer time. This is a fixed time and it would contribute only with 1/16 if level 6 
would be employed too. 


Table 1 

1 -FMG-V(1,1), 2-D, 16-EV, 5-Lev 
Residuals at Start and End of ML Cycles 



level 1 

level 2 

level 3 

level 4 

level 5 

1 

B 9 

.48E-I-2 .37E-13 

.69E+0 .97E-13 

.22E-H0 .30E-12 

.60E-1 .64E-04 

.15E-1 .14E-04 

.53E-H2 .44E-13 
.6lE-i-2 .38E-13 
.66E-I-2 .68E-13 
.56E-i-2 .39E-13 

.30E+2 .14E-12 
.30E+2 .80E-13 
.30E+2 .17E-12 
.30E-I-2 .24E-12 

.llE+2 35E-12 
.llE+2 .29E-12 
.llE+2 .30E-12 
.llE+2 .45E-12 

,30E”f' 1 .86Ey-03 
.30E+1 ,86E^03 
.30E+1 .54E-02 
.30E+1 .54E-02 

.76E-1-0 .73E-04 
.76E+0 .73E-04 
.76E+0 .llE-02 
.76E+0 .llE-02 

.52E-I-2 .12E-12 
.59E+2 .31E-12 
.61E+2 .20E-12 
.62E+2 .17E-12 
.73E+2 .13E-12 
.58E-I-2 .34E-12 

.llE+3 .32E-12 
.45E+2 .54&11 
.45E+2 .57E-11 
.45E+2 .71E-11 

.16E+2 .35E-12 
.16E+2 .32E-12 
.X6E+2 .41E-12 
.16E+2 .33Ei-12 

.44E+1 .42E-02 
.44E+1 .42E-02 
.44E+1 .39E-02 
.44E+1 .16E-02 

.llE+1 .82E:-03 
.llE+1 .82E-03 
.llE+1 .83E-03 
.llE+1 .93E-03 

•45E+2 .15E-09 
.nE+3 .41E-09 
.12E+3 .81E-11 
.12E+3 .50E-04 

.12E+3 .30E-09 
.12E+3 .16E-09 
.12E+3 .26E-09 
.12E+3 .50E-09 

.43E+2 .72E-08 
.43E-I-2 .20E-08 
.43E+2 .21E-05 
.43E+2 .16E-05 

.12E+2 .33E-02 
.12E+2 .33E-02 
.12E+2 .29E-01 
.12E+2 .29E-01 

.54E"f"2 .39E-12 
.51E+2 .70E-12 
•44E+2 .96E-12 
.53E-H2 .16E-12 


•12E‘F3 .17E-06 

.44E+2 .34E-02 




.69E+2 .19E-06 



Table 2 

1 -FMG-V(1,1), 2-D, 16-EV, 5-LEV 
Eigenvalues at End of ML Cycles 

E 1 

level 1 

level 2 

level 3 

level 4 

level 5 


.496347395806E-I-1 


.495552134150E+1 

.49550931 7773E+1 

.4954981 73425E+1 

WiTiVi 


.103677004418E-I-3 
.103677004418E+3 
.10371 E+3 
.10375 E+3 

.104634633842E+3 
.104634633842E+3 
.1046 E+3 
.1047 E+3 

.1 0487469501 2E+3 
.10487469501 2E+3 
.10491 E+3 
.10495 E+3 

.1670 E+3 

.167113893828E+3 

.167113893828E+3 

.167113893828E-I-3 

.16711389.3828E-I-3 

.16715 E+3 

.194919376181E-I-3 

.194919376181E+3 

.194962161804E+3 

.194962161804E+3 

.2024351 53808E+3 
.2024.351 53808E+3 
.202479632146E+3 
.202479632146E+3 

.204351 758395E+3 
.204351 758395E+3 
.204395 E+3 
.204396 E+3 

.204831 900326E+3 
.204831900326E+3 
.20487691 8643E+3 
.204876918643E+3 

.3291 85001 547E-I-3 
.329185001547E-I-3 
.329227787655E+3 
.329227787656E-I-3 

.38481 2002762E+3 
.38481 2002762E+3 
.3848590736 E+3 
.3848590739 E+3 

.399841022256E+3 
■399841022256E+3 
.399888846 E+3 
.399888846 E+3 

.403673103803E+3 
.4036731 03808E+3 
.403720980600E+3 
.403720980600E+3 

.248170840742E+3 

.248170840742E+3 

.248207366784E+3 

•248207366784E+3 

.424191908011E+3 

.424295844705E+3 

.483580557031E+3 

.499567983067E+3 




.32926431 3697E+3 
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Example II: Fine level Separation 

In the next example the potential V{x,y) = 5 + 35m(10x) + 2cos{l0y) causes 
a further spliting of the eigenvalues. The clusters were treated sequencially and the 
projection for each cluster was performed on the finest level to provide accurate finest 
level separation inside clusters. The results, for 9 eigenvectors, are presented in Tables 
3 and 4. A 10-FMG-V(1,1) algorithm weis used to show the constant convergence rate 
per cycle, The coarsest relaxation level for clusters 2 and 3 was level 2 and for the 
first eigenvector was level 1. On levels 1 and 2 the adaptive algorithm and few cycles 
were used. All eigenvectors came out accurately orthogonal (10 scalar products 
on level 4). 

This shows that it is enough to perform separation only within clusters. 


Table 3 

10 -FMG-V(1,1), 2>D, 9-EV, 4-Lev 
Residuals at Start and End of ML Cycles 

Evect. 

level 1 

level 2 

level 3 

level 4 

1 

2 

3 

4 

5 

6 

7 

8 
9 

0.48E+02 0.14E-13 

0.46E-I-02 0.83E-09 

0.52E-I-02 0.29E-09 

0.56E-I-02 0.56E-10 

0.54E+02 0.85E-09 
0.53E-I-02 0.57E-02 

0.53E+02 0,68E-11 
0.41E+02 0.13E-10 

0.43E+02 0.80E-02 

0.83E+00 0.12E-12 

0.30E+02 0.48E-09 

0.30E-I-02 0.12E-08 
0.30E-1-02 0.73E-09 
0.30E+02 0.55E-08 
O.llE-1-03 0.40E-05 
0.45E-h02 0.57E-05 

0.45E4-02 0.29E-05 

O.llE-1-03 0.14E-05 

0.27E+00 O.llE-11 
O.llE-l-02 0.42E-12 
O.llE+02 0.75E-12 
O.llE+02 0.93E-12 
0,llE+02 0.17E-11 
0.16E+02 0.61E-12 
0,16E-1-02 O.IOE-II 
0.16E+02 0,82E-12 
0.16E-I-02 0.83E-12 

0.72E-01 0.41E-11 

0.30E+01 0.21E-11 
0.30E+01 0.52E-11 

0.30E+01 0.56E-U 

0.30E+01 0.12E-10 

0.44E+01 0.16E-10 

0.44E+01 0.39E-10 

0.44E+01 0.33E-10 

0.44E+01 0.48E-10 




Table 4 




10 -FMG-V(1,1), 2-D, 9-EV, 4-LEV 




Eigenvalues at End of ML Cycles 


E 

level 1 

level 2 

level 3 

level 4 

1 

-0.494698319454E+01 

-0.49378951 8604E-|-01 

-0.493543833853E-r01 

-0.493481 21 4576E+01 

2 

-0.860202443918E-f02 

-0.9991 E-f-02 

-0.10367E+03 

-0.10463E+03 

3 

-0.860202443918E+02 

-0.999361E-I-02 

-0.103693lE-t-03 

-0.104650E+03 

4 

-0.860406326305E+02 

-0.999362E4-02 

-0.1036937E-I-03 

-0.104651E+03 

5 

-0.860406326305E-I-02 

-0.9997E+02 

-0.1037E-I-03 

-0.10469E+03 

6 

-0.1670E+03 

-0.19491 E-f-03 

-0.20243E-I-03 

-0.20434E+03 

7 

-0.167113893828E+03 

-0.19493E+03 

-0.20245+03 

-0.20436E+03 

8 

-0.167113893828E+03 

-0.19495E+03 

-0.20247E+03 

-0.20439E+03 

9 

-0.16713E+03 

-0.19497E-f-03 

-0.20249E+03 

-0.2044E+03 


Example III: Coarse Level Separation, 2D and 3D 

In the next two runs (Tables 5, 6) we show that a coarsest level separation {h = 
li = 1), even for clusters containing very close and degenerate eigenvalues, can be 
enough to provide accurate finest level separation. In the degenerate clusters the 
eigenvectors were not orthogonalized on finest levels but resulted so from the FMG, 
where orthogonality was imposed on the coarsest level solution (in the FMG, not 
during the cycles). This implies an 0(</iV) algorithm even for close clustered cases. 

Table 5 shows results for a problem in 2-D with a potential V{x,y) = 2 + 
0.1sm(10a: + lOy), which produces a splitting of the first cluster of four eigenval- 
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ues into two degenerate clusters whose eigenvalues are close to within 10“'’. A second 
order approximation was obtained by an 1-FMG-V(1,1) algorithm and the asymp- 
totic convergence rate per fine level cycle was 1/10. Observe the 13 equal digits of the 
degenerate eigenvalues, on all levels. On level 5, 8 cycles were performed to show the 
constant convergence rate per cycle, (see cycles 3 and 8 where the convergence rate 
is accurately 1/10). The eigenvectors came out accurately orthogonal, even in the 
degenerate eigenspaces, although the projection was performed only on the coarsest 
level, (the eigenvector’s scalar products being of order 10 on level 5). 

The same eflSciency is obtained for problems in 3-D as can be seen in Table 6. 
The potential V{x, y,z) = 2 + sm{20x d- lOt/ - 10^) determines a cluster of six close 
eigenvalues grouped into two clusters of two and four equal eigenvalues (with 13 
digits). Observe the 1/10 convergence rate in cycle 3 and the first 6 common digits 
of the eigenvalues in the cluster. 

These examples indicate that the reduction of complexity of the newly developed 
algorithm from 0{q^N) to 0{qN) is quite general and also works for rather difficult 

cases. 



1 

Table 5 || 

-FMG-V(U), 2-D, 5-EV, 5-LEV || 

iiBBai 

vectoiH 

first res T 

last res | 

eigenvalue 

Ll 

E V E L 4 II 

1 

1 

0.18E-02 

0.13E-03 

-0.19999752449715E+01 


2 

0.30E+01 

0.43E-02 

-0.101 62979203934E+03 


3 

0.30E-I-01 

0.43E-02 

-0.101 62979203934E+03 


4 

0.30E+01 

0.43E-02 

-0.101 72931140738E+03 


5 

0.30E+01 

0.43E-02 

-0.101 72931 1 40738E-f 03 


L 

E V E L 5 II 

1 


0.46E-03 

0.36E-04 

-0.19999750026388E+01 



0.76E-1-00 

0.40E-03 

-0.10186970728937E-(-03 


mM 

0.76E+00 

0.40E-03 

-0.10186970728937E-H03 



0.76E+00 

0.40E-03 

-0.10196970729590E-1-03 


5 

0.76E-1-00 

0.40E-03 

-0.10196970729590E-f03 

3 


0.35E-05 

0.33E-06 

-0.19999749801 202E+01 



0.27E-04 

0.26E-05 

-0.10186970049930E+03 


mm 

0.27E-04 

0.26E^05 

-0.101 86970049930E-H03 



0.27E-04 

0.26E-05 

-0.101 96970049780E+03 


5 

0.27E-04 

0.26E-05 

-0 . 1 01 96970049780E+03 


1 

0.97E-11 

0.26B1-10 

-0.199997497991 42E+01 



0.29E-09 

0.31E-10 

-0.101 86970048459E-h03 



0,29E-09 

0.31E-10 

-0.10186970048459E+03 



0.29E-09 

0.31E-10 

-0.10196970048302E-f03 


1 ^ 

0.29E1-09 

0.31E-10 

-0.10196970048302E-(-03 
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L_ 

1 

-FMG-V(1 

Table 6 

1),3-D, T-EV, 3-LEV, 

1 cycle II vector 

first res 

last res 

eigenvalue | 

1 LEVELS 0 

3 

1 

0.15E-02 

O.lOE-03 

-0.19991341655960E+01 


2 

0.25E-01 

0.25E-02 

-0.10072012662990E+03 


3 

0.25E-01 

0.25E-02 

-0.10072012662990E+03 


4 

0.23E-01 

0.21 E-02 

-0.1 00720682691 98E+03 


5 

0.23E-01 

0.21E-02 

-0.1 00720682691 98E-L03 


6 

0.23E-01 

0.21E-02 

-0.1 00720682691 98E-L03 


7 

0.23E-01 

0.21 E-02 

-0.1 00720682691 98E-L03 


5 Conclusions 

A robust and efficient ML algorithm to compute a few eigenvectors and the corre- 
sponding eigenvalues for large scale eigenvalue problems ha^ been developed. The 
algorithm’s robustness results from the adaptive completion and treatment of clus- 
ters, the simultaneous treatment of solutions in each cluster, and from tests which 
monitor the algorithm’s convergence and efficiency. The alg orit hm treats central dif- 
ficulties such as: the poor solution representation, on coarse levels, the existence of 
clustered eigenvalues, the approximation of incomplete clusters, and the mixing of 
approximated eigenvectors during the solution process. Its eigenvector separation ef- 
ficiency stems from a new ML projection technique which is a generalization of the 
Rayleigh Ritz projection, combined with backrotations."^ 

In the cases when the algorithm properly separates the eigenvectors on coarse 
levels, the algorithm’s complexity is of 0{qN) for q eigenvectors of size N on the 
finest level. The numerical tests showed that an accurate fine level separation was 
obtained by the coarse level projection, even for problems with very close or equal 
eigenvalues. 

The results of the numerical tests for Schrodinger eigenvalue problems, show that 
the algorithm achieved the same accuracy, using the same amount of work (per eigen- 
vector), as the Poisson multigrid solver. A second order approximation is obtained 
using the 5-point in 2-D and 7-point in 3-D discretized Laplaceian, by 1-FMG-V(1,1) 
in 0{qN) work. This means that the work was of order a few (2-4) fine level re- 
laxations per eigenvector. The adaptive work was only part of the fine level work 
and enhanced the fine level cycle efficiency. Constant convergence rate per cycle was 
obtained for the presented cases. The robustness of the algorithm has been demon- 
strated on problems with eigenvalue distribution that present special difficulties. 
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